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Microdroplet oscillations during optical pulling 

Simen A. Ellingsen 1 

Department of Energy and Process Engineering, Norwegian University of Science and Technology, N-7491 Trondheim, 
Norwa^]\ 

It was recently shown theoretically that it is possible to pull a spherical dielectric body towards the source of 
a laser beam [Nature Photonics 5, 531 (2011)], a result with immediate consequences to optical manipulation 
of small droplets. Optical pulling can be realised e.g. using a diffraction free Bessel beam, and is expected 
to be of great importance in manipulation of microscopic droplets in micro- and nanofluidics. Compared to 
conventional optical pushing, however, the radio of optical net force to stress acting on a droplet is much 
smaller, increasing the importance of oscillations. We describe the time-dependent surface deformations of a 
water microdroplet under optical pulling to linear order in the deformation. Shape oscillations have a lifetime 
in the order of microseconds for droplet radii of a few micrometers. The force density acting on the initially 
spherical droplet is strongly peaked near the poles on the beam axis, causing the deformations to take the 
form of jet-like protrusions. 



I. INTRODUCTION 

The use of optical forces in microfluidics is a field 
of research in rapid growth. Applications are already 
several, and many more seem to be promised^. The 
use of light to manipulate microscopic flows is an at- 
tractive prospect since it is contacless and nondestruc- 
tive, and easily reconfigured 3 . This marriage of op- 
tics and fluid mechanics has been dubbed optofluidics, 
and was recently the topic of a focus issue of Nature 
Photonics 4 . Optical manipulation of microscopic fluid 
systems have been studied experimentally by several 
groups, demonstrating such phenomena as manipulation 
of microdroplets in channels 5 -, jet formations due to light 
scattering* 5 , and manipulation of surfaces^ and droplets 8 
by laser-induced Marangoni flows from surface tension 
gradients. Another beautiful example is the sorting of 
microparticles in a microfluidic environment according 
to size or dielectric constant using an optical lattice^. 

Optical manipulation of spherical particles is partic- 
ularly attractive from an experimental point of view, 
not least because the force is exactly calculable with- 
out heavyweight numerics. Optical trappin g 10 ' 11 of 
particles as well as pushing with radiation pressure 12 
have been the fundamental degrees of freedom in opti- 
cal micromanipulatio n 13 ' 14 . The motivation behind the 
present study, however, is the recent discovery that it is 
also possible to pull a spherical object against the direc- 
tion of propagation of a laser beam 1 ^, adding new possi- 
bilities of optical manipulation. Application to spherical 
microdroplets is an obvious idea, yet the deformability 
of liquid droplets adds a new level of complexity to the 
theoretical description. Optical deformation of droplets 
was studied by Zhang and Chang some time ago 1 ^, and 
was the topic of a theoretical investigation by Lai et al. 
soon after—, a work generalized by Brevik and Kluge 1 ^. 

Chen et al. demonstrated 15 that optical pulling, while 
impossible with a standard plane wave source, may be re- 
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alized using a Bessel beam 1 ^, a type of optical beam with 
the additional and attractive property that it is prop- 
agation invariant and does not spread. Its usefulness 
for micromanipulation has already been demonstrated 
experimentally 2 ^—. 

In this article we consider the "optical pulling" situa- 
tion analyzed in Ref. Hsl with a view to describing fluid 
motion of a droplet in the time after it has been sub- 
jected to an optical pulling pulse. The droplet, whose 
radius must lie in one of a few intervals a few times the 
wavelength of incident light, is assumed to remain spher- 
ical during the duration of the pulse, whereupon surface 
movement set in motion by the pulse manifests itself and 
eventally dies out due to viscosity. We wish furthermore 
to compare the pulling situation to the conventional push 
by a uniform plane wave, such as may be produced by a 
beam whose width far exceeds the extent of the droplet. 

It remains an unanswered question how large droplet 
deformation can be before the optical pulling effect is de- 
stroyed, however since it is an effect of subtle diffraction 
and interference interactions within the droplet, it is nat- 
ural to presume that pulling, at least in its simplest form, 
requires the shape to be close to spherical. We therefore 
consider the case where a droplet is illuminated by a laser 
pulse which is short compared to the hydrodynamic re- 
sponse time of the droplet. This is the case studied in the 
classical experiment by Zhang and ChangiS. The pulse 
transfers an impulse to the droplet, and, in addition to 
moving the droplet as a whole, sets surface oscillations 
in motion which gradually die out due to viscosity. Once 
oscillations have vanished, a new pulse may be transmit- 
ted, and so on. 

In the next section the hydrodynamic equations of mo- 
tion for the droplet surface are laid out, where pertur- 
bations are included to linear order. Section [IIII presents 
the optical theory necessary to calculate the optical force 
density from a (general) Bessel beam as a function of po- 
lar angle acting on the initially spherical droplet, and nu- 
merical investigation of the system is considered in llVI In 
this section we assume the Bessel beam to be zeroth or- 
der and consider the case of a water droplet. The optical 
pulling case is compared to the conventional situation of 
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FIG. 1. Sketch of the set-up considered: a droplet is illu- 
minated by a Bessel beam such that the net force is a pull 
towards the laser source. A stress acts on the droplet setting 
in motion surface oscillations. 

a net pushing force before conclusions in section [Vj Var- 
ious elements of the optical theory employed are found 
in appendices. 

II. HYDRODYNAMIC EQUATIONS OF SURFACE 
MOTION 

We consider a droplet illuminated by a short laser pulse 
of duration to- As was the case in the experiment of 
Zhang and Chang we assume to to be short compared to 
the hydrodynamical response time, so that the droplet 
can be approximated as being of spherical shape for the 
duration of the pulse. This allows the use of the well 
known theory of Mie scattering, as we will do in Section 
IIIII We assume the droplet to be incompressible, and 
may therefore omit the contribution from electrostriction 
from the start^. 

A droplet irradiated by a laser beam is subjected to an 
optical force density distributed over the surface of the 
droplet. Using the completeness of Legendre polynomi- 
als, the force per unit area over the droplet surface may 
be written 

a(Q,t) =a(Cl)[e(t) -&(t-t )}; (la) 

OO I 

a(Q)=J2 vimPr (cos 8)e im *. (lb) 

1=0 m=— i 

Here P/ 71 are the associated Legendre functions, where 
in particular P ; ° = p are Legendre polynomials, and O 
is the unit step function. Q denotes a solid angle (0, </>)• 
We will calculate the details of this force density compo- 
nents ai m in the next section, but take them to be known 
constants for now. The following derivation follows that 
of Brevik and Klugei£, in turn based on that of Lai and 
co- workers^. For further details the reader may refer to 
these references. 

We will find in Section [TTT1 that er(f2) only depends on 
polar angle 0, not <j). Thus only the m = term con- 
tributes to the sum (Jib}. Clearly, an axially symmet- 
ric force will give rise to axially symmetric surface mo- 



tion, and we may simplify our notation to axial symmetry 
henceforth. We define 

OO 

*(6) = J2<n p l( cosd )- ( 2 ) 

1=0 

Let the surface of the droplet be described by the ra- 
dius function 

OO 

r{6,t) = a + ^Thi{t)Pi(cos6) =a + h(t), (3) 

where hi(t) will be our unknown functions to be solved 
in the following, a is the radius of the droplet in its 
inititial spherical form. Due to conservation of droplet 
volume, the I = solution is trivially ho — and will 
not be considered further. To linear order in hi /a all 
the modes I > 1 automatically conserve volume since 
siri 0P;>1 (cos 0) = 0. 

We turn now to the Navier-Stokes equation governing 
the fluid movement in the droplet. Working to leading 
order in the surface perturbations, convective terms are 
neglected, leaving 

f) 1 1 

— = Vp + vV 2 v+-a(6,t)5(r-a)Y. (4) 

at g g 

Here g is (constant) fluid density and v is kinematic vis- 
cosity and a hat denotes a unit vector. Incompressibility 
implies that V 2 p = in the interior of the droplet since 
no electromagnetic forces act there. Thus 

OO £ 

K0) = E(0 Pi p i( cosd )- (5) 
1=1 

The viscous term in Eq. (0| is derived e.g. in Ref. [ID 
and its radial component may be written, to leading order 
in hi I a 

OO 

v \J 2 Vr =-2^(nhi(t)Pi(cosO); (6a) 
i=i 

/ii=4(2« 2 -I -1). (6b) 

The radial component of Eq. (|4]) in the interior now reads 
as r — > a, making use of the orthogonality of the p 
expansion 

hi + 2fnhi = -—pi- 
ag 

The boundary condition at the surface is that the in- 
ward pressure Ap s t ' minus the outward optical pressure 
a balances the hydrodynamic normal stress at r = a, 
giving p\ r=a ~ Ap s t - — a. The droplet's surface tension 
creates a pressure discontinuity across the surface whose 
details were worked out by Lai et alJ^ to leading order 
in hi J a: 

Ap 8 -*' = -L + -L V(P + l - 2)hi{t)Pi{cosO) (7) 

17 I7 Z ^— ' 
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with 7 the surface tension coefficient. 
Defining 

ga 6 

we obtain the equations of motion for hi m (t): 

I 



(8) 



hi + luihi + Wj h t 



ga 



-ai[Q(t)-Q(t-t )}. (9) 



This equation has the form of a forced harmonic oscillator 
with damping coefficient pi supplied by viscosity. 

Eq © can be solved quite simply using Laplace trans- 
formation, letting hi(s) = £{hi(t)}. The Laplace trans- 
formed equation reads 

- ^-*o s 

hi(s) 



la, 



ga s(s 2 + 2fiis + u)f) ' 
whereupon inverse transformation gives the solution 



h(t) 



gauif 



| 1 - (— sin + cos 7; ^e~ rt * 6(t) 



1 - 



x e 
where 



( — sin 7 /(t - t ) + cos7/(i - to) ) 

e(t-to)} 



-jU;(t-to) 



(10) 
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Oscillations are underdamped when > pi, hence 7; 
real, which is what we consider. The oscillation ampli- 
tude for mode I is proportional to the force coefficient 
oi for that mode, as one would expect. Mode I has a 
lifetime 

- -i- 

n ~^ 1 ~u(2p-i-iy 

For water lifetimes are in the order of a microsecond for 
droplets of a few micrometers' radius. Moreover we see 
that since r ~ l~ 2 for large I, the higher order modes 
will damp out quickly, so the sum may be truncated at 
lower I as time passes, which is useful. Finally we notice 
that when I — 1 we have uj\ = p\ = 0, so Eq. ([9]) gives 
hi(t) — for t > to. The linear behaviour of h\(t) that 
this implies, and the fact that P\ (cos 6) = cos#, mean 
that the I = 1 equation describes the translative motion 
of the droplet as a whole due to pulling or pushing, not 
a surface perturbation. For our present purposes we will 
not consider this term, therefore. 

A note on electrostriction is warranted at this time. 
For an incompressible fluid the speed of sound is infinite, 
so that an electrostrictive force density due to the pres- 
ence of the laser beam is immediately compensated by a 
mechanical pressure. Electrostriction therefore plays no 
role when it comes to the actual motion of the surface. It 
is of importance, however, to the droplet's stability: the 
inwardly directed electrostrictive force is larger than the 
outwardly directed net optical force. For further discus- 

c.;™ 23.24 



III. OPTICAL FORCE DENSITY 

With the equation of motion solution of Eq. (fTU|) . what 
is required in order to determine the droplet surface os- 
cillations are the coefficients <t; of the surface force den- 
sity expansion JTJ). The optical force density is obtained 
from Mie scattering theory 2 ^. We make use of the no- 
tation of Barton et ali^£, to whose exposition the reader 
should refer for further details. In the following we de- 
note spherical coordinates (r, 9, (j>) with 9 the polar angle 
and origin at the droplet centre, while cylindrical coordi- 
nates are denoted (p, <p, z). As numerical benchmarks we 
have reproduced the force calculations of Refs.ll5l andl27l 
for the Bessel beam and plane wave, respectively. 

We consider the case of an incident field S l imping- 
ing on an initially spherical droplet with refractive index 
ni, embedded in a medium of refractive index U2- In 
the present paper we shall assume ri\ and ni to be real 
for simplicity, although the generalisation to a droplet of 
complex m, hence absorbing droplet, is straightforward. 
Note that the realisation of optical pulling depends upon 
the imaginary part of n± being small since the droplet 
receives a pushing impulse for every impinging photon 
whose momentum it absorbs, eventually destroying the 
pulling effect. 

The radial force density acting on a surface element of 
the sphere may be calculated by integrating the gradient 
of the Maxwell stress tensor 

T = S®T) + n®B-\{S-T> + n-B)l (11) 

across the surface, 

cr(Q) = (a rr ) = (T„(r = a+) - T rr (r = a")). (12) 

Here 1 is the unit matrix and (• • • ) denotes time aver- 
age over an optical period. (For an isotropic, dielectric 
medium the tensors of Minkowski and Abraham coincide. 
The Abraham forced oscillates out and gives no contri- 
bution in the optical case.) We let S = Rc{Ee 4W *}, etc., 
where E, D, B and H are complex field vectors. For field 
quantities Xi and Xj, we have (XiXj) = jRcjXjX*}, 
and in particular (Xf) — ^\Xi\ 2 . 

The simplest option now is to express the stress tensor 
both outside and inside the spherical surface in terms 
of the interior fields, since the external fields have both 
incident and scattered components. Using the continuity 
of D ri Eg,Ecf, and H across the surface, the force density 
may be written 



a(Q) 



e «2 



(n'-l)(n'\E?\' + \E?\' + \E%n (13) 



Superscript w signifies that the fields are evaluated just 
within the surface of the sphere, at radius r = a~ , and 
we have defined the shorthand 



ni 

712 



(14) 



4 




FIG. 2. Polar plot of the angular distribution of force density a plotted as 5 + a(8)/a where a = \ f£ d8 sin8a(8). The shading 
indicates the value of a (9) /a. Comparison of pulling with order zero Bessel beam (net force towards the left) and conventional 
pushing by plane wave (net force towards the right), a = 11.35 in both cases. 



The internal fields are related to the incident electric 
and magnetic fields from the Bessel beam via the rela- 
tions quoted in Appendix [B] For further details, refer 
to Refs. [lU and [26|. The incident electromagnetic field 
components of an order m Bessel beam are given in Ap- 
pendix |XJ For comparison we consider a circularly po- 
larised plane wave, the formalism for which is quoted in 
Appendix [Cj The force density as it is distributed as a 
function of polar angle 9 is plotted in figure [21 to which 
we will return for discussions later. 

A Bessel beam's plane wave components form cones 
of angle #0 with the direction of propagation, so that 
k z = fccos#o and k± = fcsinfo- A factor e~ tult is omitted 
in all field components here and henceforth. 

Internal field quantities for insertion into Eq. (|13l) are 
calculated from the incident field quantities as laid out 
in appendix |Bj We require the quantities Ai m and B\ m 
defined in Eq. (|B3[) . which we find as 



Aim 
Blm 



1 

n 2 



V4tt(2/ + 1) / (I -to)! A 
1(1 + 1) V (l + m)\ Bi m 



! in 



where we have defined 



(A,B) lm = 
cos a z u 



4>i(a) 



sina,M 



± 



smci,M 



cos a z u 
ima 



— J m (b±)r]TM,TE ■ , 



Vl - U 2 \pr(u) for I + to 



even 
odd 



J m {b±)uriTM,TE 

Jm(b±)r]TE,TM 

(15) 



(before and after comma pertains to Ai m and Bi m , re- 
spectively; upper and lower pertains to even/odd values 
of I + to). Herein, k — 2tt/X — ri2o;/c, a Z! j_ = fc Zi j_a, 
2 and 



b± 



.vT 



= ka = 



2na 



(16) 



Moreover J' m (x) — -§^Jm(x) where J m is the cylindrical 
Bessel function of the first kind. 7/te and ryxM are the 
relative weights of TE (transverse electric, w.r.t. direc- 
tion of propagation) and TM (transverse magnetic) po- 
larizations of the incident beami 5 -. These weights in gen- 
eral carry a relative phase S between the polarizations 15 , 
6 = Arg (?7tm / ?7te ) + "7 2 - 

We make the convenient definitions 



ci m =iAi m [nipi(na)^' (a) - ^'(ria)( ( (1) (a)] 1 (17a) 
dim =iB im [V«(™*)£, (1) '(a) - n^na)^ (a)}- 1 (17b) 



wherewith, using the Mie theory in appendix [B] we may 
finally write down the electromagnetic force density on 
the droplet surface from an order toq Bessel beam as 



fa\ ^/pa 2 ,. 2 n (2fc + l)(2l + l)(k - m Q )\(l - to )! f k(k + 1)1(1 + 1) . 
^ )= ^- (n ' 1} )^ i k(k + mi + l)(k + mo )Kl + mo)\ \ a* c kmo Ci mo Mna)Mna)P k P t 

+ [ckm CimJ'k(naWi(na) + d km Jl m J k (na)ipi(na)](siv? eP^'P^ ' + ml-± 



pmo pmo 

hi- 1 rn ? 
k 1 I + TO . 2 a 

sin 



+ m [c kmo d^ mo ip' k (na)ipi(na) + d kmo 4 mo ^k(na)ij[(na)][P 7 k na P/™ ]'} (18) 



where we suppress the argument cos# of Legendre func- tions, P, m/ (a;) = dP; m (a;)/da;, and * denotes complex con- 



jugate. Iq is the central intensity in the case mo = 0. 

As promised in Section UH a only depends on the po- 
lar angle 9 since the incident field components depend on 
<p only through an overall phase factor exp(iTO O 0) which 
vanishes in Eq. (I13[) . This means that only m = con- 
tributes to the expansion in Eq. (p}, which reduces to 
Eq. ©. 

From definition <j2j) and orthogonality relation (|B6|) it 
follows that the coefficients 07 are found as 

07 = (Z + ±)/ d9shi9a(0)Pi{cos9). (19) 
Jo 

Formally, ai requires calculation of a doubly infinite sum, 
yet most of the integrals resulting from inserting (|18[) into 
(Q33) are zero. 

IV. NUMERICAL INVESTIGATION 

For numerical purposes, let us concentrate on the sim- 
plest case of a TM polarized Bessel beam of order m Q = 0. 
This is the case considered in Fig. 1 of Ref. HH, and im- 
plies ?7tm = l,r?TE = 0. We let 9 = 78.5°. The far 
field scattering off a sphere from an order zero beam was 
recently calculated by Mitri2&. 

Numerically one finds that the highest value of I that 
must be included in expansion §5§ is of order a, as is 
physically reasonable since the wavelength of the incom- 
ing light is the smallest spatial scale of the surface force, 
and a is the number of wavelengths in a circumference. 
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FIG. 3. Nondimensionalised optical force per cross-section 
area, Eq. (|20[) , on front half and rear half of a water droplet 
during (a) conventional pushing by plane wave and (b) push- 
ing/pulling by TM order zero Bessel beam. The middle graph 
in both panels shows the net force on the sphere. We let 
A = 1.064/^m. 



A. Net force versus stress on the droplet 

Before considering the oscillations resulting from the 
optical pulling force, we will regard the integrated force 
acting on the front and back halves of the droplet and see 
how they compare to the same situation for a circularly 
polarized incident plane wave. This gives an intuitive 
picture of how the net propulsion force (obtained when 
integrating over the full sphere) compares with the opti- 
cal stress tending to pull the droplet apart. Field equa- 
tions for such a plane wave are found in Appendix [Cj The 
force in z direction acting on the front and back halves of 
the sphere are found by integrating over the appropriate 
solid angles. Working in terms of the dimensionless force 



Q = 



c(F z 



na 2 n2 



'-■ (F z )=a 2 / dna(9) 



(20) 



half-sphere forces are (mo = 0, TM polarization) 

oo 

Qi =(n 2 - 1) 2J c k Ci[ip k (na)ipi{na)V t 



k,l=l 



(21) 



with i = {front, back}, and 



V^ nt =(2fc + l)(2i + 1) J d«P,(u)F ; ( M ); 

and corresponding quantities for the back half by replac- 
ing the integrals with f du. 

Figure [3] shows how the z-directed integrated surface 
forces acting on either half of the sphere compare with 
the net pulling or pushing force. Clearly a liquid droplet 
does not behave like two rigid hemispheres, but the graph 
nevertheless gives an intuitive picture of the magnitude 
of optical stress pulling the droplet apart as compared 
to the propulsion force. It is immediately clear that in 
regions where optical pulling is possible, i.e., where the 
net force in the Bessel beam case is negative, the stress 
is an order of magnitude larger than the net force. This 
is in contrast with the plane wave case where the two are 
in the same order of magnitude. Thus it is clear that in 
order to obtain a pulling force as large as the conventional 
pushing force, droplet oscillations set in motion by the 
optical stress will be comparatively greater. [Note that 
the difference in the numerical values of Q in the two 
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cases in Fig. [3] is of little importance; to obtain the same 
irradiated power, field amplitude E must be larger for 
the Bessel beam than the plane wave, accounting for the 
difference! 



B. Droplet oscillations 

For numerical purposes we take the values for water— 
and let the surrounding medium be air (712 = 1) 



g =997kg/ni 
rii =1.33; 



7 = 0.073N/m; 
v = 1.01/iN/m. 



As in Rcf. 1 — we choose A = 1.064/xm. For optical pulling 
we choose the fourth pulling regime from the left in figure 
[3}d, to wit 

a = 1.922/wn; a = 11.35. 

The reason for not choosing the smallest a value is to 
obtain a more interesting oscillation structure which is 
not present when a is of order unity or below. Note how- 
ever that the stress relative to net force is much greater 
for the smallest optical pulling radius interval (around 
a = 0.5fim in figure (3})), so that oscillations will be ex- 
pected to be more violent in that case, albeit somewhat 
compensated by the greater influence of surface tension 
for a smaller droplet. 

In order to compare the Bessel beam pulling case with 
the conventional plane wave case, we let the two cases 
have the same laser power, i.e., the Poynting vector inte- 
grated over area it a 2 should be the same for the two cases. 
For the plane wave, W w = E^Ae ikz with A = x + iy. 
The intensities are then related by 



1 o 



cot 9 



. J ( a -i-) Jo(a±)Ji(a!±) + Ji(a±) . 

1q sin Wo L a± J 

(22) 

With the numbers used in the examples, Iq is then greater 
than Jg W by a factor 88.2. 

For the Bessel beam we require numerical calculation 
of the coefficients [these are also analytically calculable, 
but a numerical treatment was found to be simpler] 



5/ m „=(^ + l)(2m + l)(2n + l) / AuJ^u) J m {u) J n {u); 

(23a) 

2 _ (2/ + l)(2m + l)(2n+l) 
lmn ~ m(m + l)n(n + l) 

x j du{\ - u 2 )Ji( u )J^(«)j;(«). (23b) 

Since we consider only moderate values of a, numerical 
calulation of all coefficients could be done within a rea- 
sonable time. 

We deliberately consider a rather high-intensity laser 
beam so as to visualize the oscillations most clearly. As 



intensity in the plane wave case we use Iq W = 2W//xm 2 
(i.e., a laser power of ), which is about a fifth of that used 
in the Zhang and Chang experiment and of that as- 
sumed by Brevik and Klugc 1 ^ as well as Lai et alJi. The 
resulting oscillations are depicted in Fig. [4] for different 
times after the onset of the laser pulse, which we assume 
to have duration to = 0.04^s and 0.02/is. The two differ- 
ent durations show the effect of reducing or lenghening 
pulse time. The oscillations are then somewhat larger 
than might reasonably be described by the linear the- 
ory herein, but nevertheless form an instructive example 
in that oscillations are clearly visible. The oscillatory 
part of the radius function ([3]) is directly proportional to 
la, hence the shapes of the perturbations as a function 
of time are the same also for smaller intensities where 
linear theory is accurate. The droplet has returned to 
approximately spherical shape after 1-2/^s. Clearly, any 
adverse effects of oscillations can be reduced by lowering 
laser intensity, at the cost of a smaller pulling force. 

We find that, especially in the case of to — 0.04/is, the 
droplet is deformed significantly during the duration of 
the pulse with our numbers. It is not known at what 
point deformations destroy the optical pulling effect, yet 
it seems likely that at this intensity level a shorter pulse 
might be required in practice. A similar net momen- 
tum may transferred to the droplet if one could use a 
shorter pulse but at higher intensity. In Fig. [5] this is 
shown. Here, surface deformation is smaller at the the 
time of pulse turnoff, while the net momentum trans- 
ferred (which is proportional to io^o) is similar. 

In order to compare the case of optical pulling with 
conventional pushing we calculate the corresponding 
droplet oscillations when the droplet is subjected to a 
pulse of the same length and power, but now from a cir- 
cularly polarized plane wave beam of infinite width. The 
formalism of such plane waves is laid out in Appendix [C] 
(c.f. also Refs. GJ, [U and f30h . The result is shown in 
figure M for I$ w = 2W//xm 2 and t = 0.04^s. Although 
the net pushing force in this case is more than 6 times 
greater than the pulling force in the corresponding case 
in figure |4l droplet deformation is much smaller. 

The reason for this is two-fold: firstly and most im- 
portantly it is due to the much greater stress-to-net force 
ratio in this shown in Fig.[3p, but a second effect 

is a relatively higher concentration of force density near 
the symmetry axis, which makes the initial perturbations 
take the form of jet-like protrusions near the droplet's 
poles. One may see this quite clearly in the plot of the 
force distribution in Fig. [2] It is possible that the "jet" 
protrusion deformation effect can be reduced by using 
a higher order Bessel beam or a combination of Bessel 
modes; a zeroth order Bessel beam has its peak inten- 
sity near the beam center whereas higher order beams 
show peaks at non-zero p. We plan to investigate such 
possibilities in a future study. 
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FIG. 4. Droplet cross sections during optical pulling by a m = Bessel beam, pulse duration 0.04/xs (solid line) and 0.02/is 
(dashed line). The laser beam enters from the left. See main text for details. 
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FIG. 5. Same as figure [4] but with shorter pulse duration (to = 0.005^is) and higher intensity (J pw = lOW/^tm 2 ). 



V. CONCLUSIONS 



We have investigated the surface oscillations on a 
droplet under optical pulling by a zeroth order Bessel 
beam and made comparison to the case of pushing by 
a simple plane wave. The possibility to pull objects to- 
wards the source of a laser beam was predicted recently^, 
and is expected to be of great technological interest in 
microfluidics. Compared to conventional pushing of a 
droplet, the stress acting on the droplet is much greater 
under pulling, setting surface oscillations in motion. 
Thus even though the droplets in question are small (a 
few micrometers' radius), surface deformations may have 
to be taken into account even at laser intensities where 
a pushing pulse results in no noticeable deformation. 



Other than the fact that droplet stress-to-net force ra- 
tio is much higher during optical pulling than pushing, 
additional enhancement of deformations stem from the 
relatively higher concentration of force density near the 
beam center. The force density acting on the initally 
spherical droplet is concentrated near the poles at 8 = 
and 9 = w, causing jet-like protrusions to appear along 
the symmetry axis. 



Relaxation times for the oscillations due to viscosity 
are furthermore discussed. For a water droplet of ra- 
dius in the order of a few micrometers, the droplet has 
returned to spherical shape after 1-2/zs. 
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f=0.1|js 




0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 

Time from start of pulse (fjs) 

FIG. 6. Surface oscillations irradiated by same power as in 
figure 21 but with a circularly polarized plane wave. Panel 

(a) shows approximately the maximum perturbation, panel 

(b) shows the relative elevation h(t)/a as a function of time 
at back (9 — 0) and front (0 — n) of droplet. 
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The magnetic field components are found from Maxwell's 
equations as 

h\= J mo (k±P) + ; J mo {k±P); (A3a) 



hi 



p 

ifc z TO /3 T E 
kp 



J ma {k±p) + 'tf™k±J' mo {kxp). (A3b) 

(A3c) 



7 J mo (k±p). 



Here (3te — ifc?7TE/&j_, Ptm — kr/TM/k^. Note that 
Eo is not the complex amplitude, but is chosen so that 
(£ 2 ) = Eq at p — for the case mo = 0, TM polarization. 
Central intensity is then Iq = eq^cEq. 
In spherical coordinates one readily finds 



ee 



TE 



ifc z fcj_/3 TM , \ sin 
Jmo+ k mo J cos 9 



, /?TMfcj 



COS0 

" t0 sin 



(A4a) 



and 



'toq/3tm 



■J,, 



k z k ± (3- 



TE 



J m 



sin (9 
cos 6* 



=F 



i/^TE^l 



COS# 



(A5a) 



Note furthermore that z = rcosO in Eq (j Al[) and p = 
rsin#. Component E\ is naturally the same in both 
systems. 



Appendix B: Mie theory relating incident and internal fields 

The internal electric field components are expanded in 
spherical harmonics according tc 26 



oo / 



Appendix A: Incident and internal fields from a Bessel beam 

The incident electric field £ from a Bessel beam of 
order m-o is given in cylindrical coordinates (p, 4>, z) asi^ 



H 1 n 2 Ho iv v ' 



with Hq = Eo/(cpo) and 



e i=^UM + ^4(Mi (A2a) 



k z m l3 T M 
kp 



J mo (k±p) - pTEk±J'(k ±P ); (A2b) 



^ = — ; J mo (k±p). 



(A2c) 



£™=£ ^ ^ I(I + l)Q m Vi(na)3Mn) ( Bla ) 



Z=l m=—l 

oo i 



E% =aE Y, [nci m ^{na)d d Y lm {n) 



1=1 m=—l 
dim , ,- sYl m (£l) 

mipi{na) — — — 

ri2 sin 



E™ =iaE o y^ V* \rnnci m i>[(na)-^ 
v t-^ L si 



/— 1 m=— 2 



sin 



"2 



■-0; (na)9eF; m (0) 



where = 



and the coefficients 



(Bib) 



c; m =iA; m [n 2 0/(nQ;)^ 1) '(Q;) - nipi(na)£i] 1 



(Blc) 



(B2a) 
(B2b) 
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and the incident field E 1 is contained in the quantities 

t'wrh^I ^ YM ^ (B3) 

where incident fields are evaluated at r = a + and the 
integral is over all solid angles. 
Note that 



«-< n >- f^$T$ Pric ° s ' ,)el "* <B4) 



and the orthogonality relation 



(B5) 



so 



2 (1 + mV 

d9 sin 9Pr (cos 9)P l T (cos 9) — ^Zgjl fo . ( B6 ) 



Appendix C: Circularly polarized plane wave 

A circularly polarized plane wave propagating along 
the z direction may be expressed as (Ref. H3, section 
10.3) 



i=i 

+ iy x ji{kr)X n ] 



(CI) 



with 



i VK* + 1 



= (r x V)F im (fi) (C2) 



where fc = ri2W / c and A = x + iy. The radial component 
may then be written (Ref. [30, section 10.4) at r = a + 



E 

K = -^Y «' + V4^a + l)(2/ + l)^(«)^i(^). (C3) 
or L — ' 

1=1 

The radial magnetic component is now found from 
Maxwell's equations as H l r = n-iE % r / '(i/xoc) . By using for- 
mulas (IB3I) and (IC3I) we find 



- _ 4tt(2Z + 1) _ i - 

^im — 7T\ , -.x "ml — &h 



1(1 + 1) 



71 2 



(C4) 



and the field components (|B1[) may be written 



z=i 



+ dii>i(na) — 

sm9 

ie^Ep ^ 21 + 1 
+ di<ipi(na)Pi sin( 



^(2/ + l) Q1 /. ; (na)P / 1 (C5a) 
ciipi(na)Pl' sin (9 

(C5b) 



£J(l + i) 



pi 



(C5c) 



having suppressed the argument cos (9 of the Legendre 
polynomials. Intensity is Iq W = eo^cEp. 
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